// In this gradient calculation routine we don't divide
// the gradient values by voxelSpacing. This shader 
// should be used with 2D transfer function, because
// in this case the lengths of the gradients in the 2D 
// transfer function histogram and in the shader are the 
// same. Also note that in this case the voxelSpacing 
// should be approximately equal to 1/(dim. of volume), 
// for instance for 256x256x256 volume it should be 
// approx. 0.004
// With OpenGL 3.0 capable hardware we will be able to 
// calculate the proper voxelSpacing automatically, but 
// for now we have to leave this to the user.

vec4 dataInterpolation(sampler3D data, const vec3 texCoordinates);

vec3 gradient( sampler3D data, const float scaleValue, const vec3 texCoord, const float voxelSpacing )
{
      vec3 delX = vec3(voxelSpacing, 0.0, 0.0);
      vec3 delY = vec3(0.0, voxelSpacing, 0.0);
      vec3 delZ = vec3(0.0, 0.0, voxelSpacing);
      //
      // Sample the texture
      //
      vec3 newVector = vec3(gl_TexCoord[0].xyz);
      vec3 newVectorX = newVector + delX;
      vec3 newVectorY = newVector + delY;
      vec3 newVectorZ = newVector + delZ;
      vec4 cx = dataInterpolation(data, newVectorX);
      vec4 cy = dataInterpolation(data, newVectorY);
      vec4 cz = dataInterpolation(data, newVectorZ);
      //
      // Calculate the gradient vector at the sample position.
      // 	
      vec3 grad;
      grad.x = (cx.a-scaleValue);
      grad.y = (cy.a-scaleValue);
      grad.z = (cz.a-scaleValue);
      return grad;
} 